library(papaja)
## Loading required package: tinylabels
library(tinylabels)
library(ds4ling)
## 
##  ds4ling loaded
##  Happy coding!
library(tidyverse)
## ── Attaching packages ─────────────────────────────────────── tidyverse 1.3.2
## ──
## ✔ ggplot2 3.4.0     ✔ purrr   1.0.1
## ✔ tibble  3.1.8     ✔ dplyr   1.1.0
## ✔ tidyr   1.3.0     ✔ stringr 1.5.0
## ✔ readr   2.1.3     ✔ forcats 1.0.0
## ── Conflicts ────────────────────────────────────────── tidyverse_conflicts() ──
## ✖ dplyr::filter() masks stats::filter()
## ✖ dplyr::lag()    masks stats::lag()
library(untidydata)
library(here)
## here() starts at /Users/alejandrojaume/Desktop/stats_paper
library(car)
## Loading required package: carData
## 
## Attaching package: 'car'
## 
## The following object is masked from 'package:dplyr':
## 
##     recode
## 
## The following object is masked from 'package:purrr':
## 
##     some
my_data <- read_csv(here("data", "data_raw_proyecto.csv"))
## New names:
## Rows: 180 Columns: 14
## ── Column specification
## ──────────────────────────────────────────────────────── Delimiter: "," chr
## (2): sex, conservatorio dbl (10): participant, c_nat, c_soc, ef_fi, leng_cas,
## leng_cat, leng_est, ma... lgl (2): ...13, ...14
## ℹ Use `spec()` to retrieve the full column specification for this data. ℹ
## Specify the column types or set `show_col_types = FALSE` to quiet this message.
## • `` -> `...13`
## • `` -> `...14`
view(my_data)

glimpse(my_data)
## Rows: 180
## Columns: 14
## $ participant   <dbl> 1, 2, 3, 4, 4, 6, 7, 8, 9, 10, 11, 12, 13, 14, 15, 16, 1…
## $ sex           <chr> "male", "male", "male", "female", "female", "male", "fem…
## $ conservatorio <chr> "no", "no", "no", "no", "no", "no", "no", "no", "no", "n…
## $ c_nat         <dbl> 9, 8, 8, 8, 3, 8, 9, 4, 5, 5, 5, 8, 9, 7, 8, 5, 7, 5, 5,…
## $ c_soc         <dbl> 9, 8, 10, 8, 3, 9, 9, 5, 7, 6, 6, 10, 10, 8, 10, 7, 8, 5…
## $ ef_fi         <dbl> 9, 9, 9, 9, 9, 9, 8, 8, 6, 9, 9, 9, 10, 9, 10, 10, 7, 7,…
## $ leng_cas      <dbl> 9, 7, 7, 7, 4, 9, 8, 5, 6, 5, 6, 7, 10, 8, 8, 6, 8, 6, 6…
## $ leng_cat      <dbl> 8, 7, 7, 6, 5, 9, 8, 5, 6, 6, 4, 8, 10, 6, 6, 5, 6, 6, 6…
## $ leng_est      <dbl> 9, 10, 9, 5, 3, 10, 10, 6, 5, 7, 7, 9, 10, 9, 9, 8, 9, 8…
## $ mat           <dbl> 10, 6, 5, 6, 4, 7, 9, 3, 3, 6, 4, 6, 10, 5, 7, 1, 5, 4, …
## $ plast         <dbl> 10, 8, 8, 6, 4, 9, 9, 6, 6, 10, 7, 8, 10, 7, 8, 7, 9, 5,…
## $ tecn          <dbl> 9, 7, 7, 7, 3, 9, 9, 5, 5, 7, 6, 7, 9, 6, 7, 5, 6, 4, 5,…
## $ ...13         <lgl> NA, NA, NA, NA, NA, NA, NA, NA, NA, NA, NA, NA, NA, NA, …
## $ ...14         <lgl> NA, NA, NA, NA, NA, NA, NA, NA, NA, NA, NA, NA, NA, NA, …
my_data_tidy <- select(my_data, participant, sex, conservatorio, c_nat, c_soc, ef_fi, leng_cas, leng_cat, leng_est, mat, plast, tecn)

my_data_tidy_proyecto <- my_data_tidy %>%
  write_csv(., path = "data/data_tidy_proyecto.csv")
## Warning: The `path` argument of `write_csv()` is deprecated as of readr 1.4.0.
## ℹ Please use the `file` argument instead.
summary(my_data_tidy)
##   participant         sex            conservatorio          c_nat       
##  Min.   :  1.00   Length:180         Length:180         Min.   : 2.000  
##  1st Qu.: 45.75   Class :character   Class :character   1st Qu.: 6.000  
##  Median : 90.50   Mode  :character   Mode  :character   Median : 7.000  
##  Mean   : 90.49                                         Mean   : 7.106  
##  3rd Qu.:135.25                                         3rd Qu.: 8.000  
##  Max.   :180.00                                         Max.   :10.000  
##      c_soc            ef_fi           leng_cas         leng_cat     
##  Min.   : 2.000   Min.   : 4.000   Min.   : 3.000   Min.   : 3.000  
##  1st Qu.: 6.000   1st Qu.: 7.000   1st Qu.: 6.000   1st Qu.: 6.000  
##  Median : 8.000   Median : 8.000   Median : 7.000   Median : 7.000  
##  Mean   : 7.439   Mean   : 8.039   Mean   : 6.911   Mean   : 6.794  
##  3rd Qu.: 9.000   3rd Qu.: 9.000   3rd Qu.: 8.000   3rd Qu.: 8.000  
##  Max.   :10.000   Max.   :10.000   Max.   :10.000   Max.   :10.000  
##     leng_est           mat             plast             tecn       
##  Min.   : 3.000   Min.   : 1.000   Min.   : 4.000   Min.   : 3.000  
##  1st Qu.: 6.000   1st Qu.: 5.000   1st Qu.: 7.000   1st Qu.: 6.000  
##  Median : 8.000   Median : 7.000   Median : 8.000   Median : 8.000  
##  Mean   : 7.544   Mean   : 6.828   Mean   : 7.739   Mean   : 7.483  
##  3rd Qu.: 9.000   3rd Qu.: 9.000   3rd Qu.: 9.000   3rd Qu.: 9.000  
##  Max.   :10.000   Max.   :10.000   Max.   :10.000   Max.   :10.000
my_data_tidy_proyecto %>%
  group_by(., conservatorio) %>%
  summarize(., mean_c_nat = mean(c_nat), sd_c_nat = sd(c_nat),
               mean_c_soc = mean(c_soc), sd_c_soc = sd(c_soc),
               mean_ef_fi = mean(ef_fi), sd_ef_fi = sd(ef_fi),
               mean_leng_cas = mean(leng_cas), sd_leng_cas = sd(leng_cas),
               mean_leng_cat = mean(leng_cat), sd_leng_cat = sd(leng_cat),
               mean_leng_est = mean(leng_est), sd_leng_est = sd(leng_est),
               mean_mat = mean(mat), sd_mat = sd(mat),
               mean_plast = mean(plast), sd_plast = sd(plast),
               mean_tecn = mean(tecn), sd_tecn = sd(tecn))
## # A tibble: 2 × 19
##   conservatorio mean_c…¹ sd_c_…² mean_…³ sd_c_…⁴ mean_…⁵ sd_ef…⁶ mean_…⁷ sd_le…⁸
##   <chr>            <dbl>   <dbl>   <dbl>   <dbl>   <dbl>   <dbl>   <dbl>   <dbl>
## 1 no                6.39    1.68    7.32    1.80    7.92    1.36    6.42    1.49
## 2 si                7.82    1.40    7.56    1.51    8.16    1.18    7.4     1.50
## # … with 10 more variables: mean_leng_cat <dbl>, sd_leng_cat <dbl>,
## #   mean_leng_est <dbl>, sd_leng_est <dbl>, mean_mat <dbl>, sd_mat <dbl>,
## #   mean_plast <dbl>, sd_plast <dbl>, mean_tecn <dbl>, sd_tecn <dbl>, and
## #   abbreviated variable names ¹​mean_c_nat, ²​sd_c_nat, ³​mean_c_soc, ⁴​sd_c_soc,
## #   ⁵​mean_ef_fi, ⁶​sd_ef_fi, ⁷​mean_leng_cas, ⁸​sd_leng_cas
knitr::kable(my_data_tidy_proyecto %>%
  group_by(., conservatorio, sex) %>%
  summarize(., mean_c_nat = mean(c_nat), sd_c_nat = sd(c_nat),
               mean_c_soc = mean(c_soc), sd_c_soc = sd(c_soc),
               mean_ef_fi = mean(ef_fi), sd_ef_fi = sd(ef_fi),
               mean_leng_cas = mean(leng_cas), sd_leng_cas = sd(leng_cas),
               mean_leng_cat = mean(leng_cat), sd_leng_cat = sd(leng_cat),
               mean_leng_est = mean(leng_est), sd_leng_est = sd(leng_est),
               mean_mat = mean(mat), sd_mat = sd(mat),
               mean_plast = mean(plast), sd_plast = sd(plast),
               mean_tecn = mean(tecn), sd_tecn = sd(tecn)))
## `summarise()` has grouped output by 'conservatorio'. You can override using the
## `.groups` argument.
conservatorio sex mean_c_nat sd_c_nat mean_c_soc sd_c_soc mean_ef_fi sd_ef_fi mean_leng_cas sd_leng_cas mean_leng_cat sd_leng_cat mean_leng_est sd_leng_est mean_mat sd_mat mean_plast sd_plast mean_tecn sd_tecn
no female 6.088889 1.635157 7.000000 1.651446 7.888889 1.569919 6.111111 1.480513 5.977778 1.514909 7.066667 1.737291 5.533333 2.051607 6.955556 1.476414 6.577778 1.514909
no male 6.688889 1.689839 7.644444 1.896834 7.955556 1.127256 6.733333 1.452271 6.444444 1.574737 7.600000 1.875682 6.244444 1.967411 7.600000 1.543314 7.222222 1.490712
si female 7.777778 1.222681 7.644444 1.539710 8.022222 1.157758 7.422222 1.469213 7.377778 1.434777 7.777778 1.636083 7.622222 1.613704 8.177778 1.211477 8.066667 1.268499
si male 7.866667 1.575379 7.466667 1.501514 8.288889 1.198905 7.377778 1.541677 7.377778 1.511906 7.733333 1.601136 7.911111 1.794211 8.222222 1.020002 8.066667 1.250455
my_data_tidy_proyecto %>%
  ggplot(aes(x = conservatorio, y = c_nat, color = sex)) + 
  geom_boxplot()      

my_data_tidy_proyecto %>%
  ggplot(aes(x = conservatorio, y = c_soc, color = sex)) + 
  geom_boxplot()

my_data_tidy_proyecto %>%
  ggplot(aes(x = conservatorio, y = ef_fi, color = sex)) + 
  geom_boxplot()

my_data_tidy_proyecto %>%
  ggplot(aes(x = conservatorio, y = leng_cas, color = sex)) + 
  geom_boxplot()

my_data_tidy_proyecto %>%
  ggplot(aes(x = conservatorio, y = leng_cat, color = sex)) + 
  geom_boxplot()

my_data_tidy_proyecto %>%
  ggplot(aes(x = conservatorio, y = leng_est, color = sex)) + 
  geom_boxplot()

my_data_tidy_proyecto %>%
  ggplot(aes(x = conservatorio, y = mat, color = sex)) + 
  geom_boxplot()

my_data_tidy_proyecto %>%
  ggplot(aes(x = conservatorio, y = plast, color = sex)) + 
  geom_boxplot()

my_data_tidy_proyecto %>%
  ggplot(aes(x = conservatorio, y = tecn, color = sex)) + 
  geom_boxplot()

mod_null_c_nat <- lm(c_nat ~ 1, data = my_data_tidy_proyecto)

mod1_c_nat <- lm(c_nat ~ 1 + conservatorio, data = my_data_tidy_proyecto)

mod2_c_nat <- lm(c_nat ~ 1 + conservatorio + sex, data = my_data_tidy_proyecto)

mod_int_c_nat <- lm(c_nat ~ 1 + conservatorio + sex + conservatorio:sex, data = my_data_tidy_proyecto)

summary(mod_int_c_nat)
## 
## Call:
## lm(formula = c_nat ~ 1 + conservatorio + sex + conservatorio:sex, 
##     data = my_data_tidy_proyecto)
## 
## Residuals:
##     Min      1Q  Median      3Q     Max 
## -4.8667 -0.8667  0.1333  1.2222  2.9111 
## 
## Coefficients:
##                         Estimate Std. Error t value Pr(>|t|)    
## (Intercept)               6.0889     0.2298  26.496  < 2e-16 ***
## conservatoriosi           1.6889     0.3250   5.197 5.58e-07 ***
## sexmale                   0.6000     0.3250   1.846   0.0665 .  
## conservatoriosi:sexmale  -0.5111     0.4596  -1.112   0.2676    
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## 
## Residual standard error: 1.542 on 176 degrees of freedom
## Multiple R-squared:  0.1941, Adjusted R-squared:  0.1803 
## F-statistic: 14.13 on 3 and 176 DF,  p-value: 2.716e-08
anova(mod_null_c_nat, mod1_c_nat, mod2_c_nat, mod_int_c_nat)
## Analysis of Variance Table
## 
## Model 1: c_nat ~ 1
## Model 2: c_nat ~ 1 + conservatorio
## Model 3: c_nat ~ 1 + conservatorio + sex
## Model 4: c_nat ~ 1 + conservatorio + sex + conservatorio:sex
##   Res.Df    RSS Df Sum of Sq       F    Pr(>F)    
## 1    179 518.99                                   
## 2    178 426.54  1    92.450 38.9015 3.216e-09 ***
## 3    177 421.21  1     5.339  2.2465    0.1357    
## 4    176 418.27  1     2.939  1.2366    0.2676    
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
plot(mod_int_c_nat)

durbinWatsonTest(mod_int_c_nat)
##  lag Autocorrelation D-W Statistic p-value
##    1      0.04504658      1.878689   0.368
##  Alternative hypothesis: rho != 0
mod_null_c_soc <- lm(c_soc ~ 1, data = my_data_tidy_proyecto)

mod1_c_soc <- lm(c_soc ~ 1 + conservatorio, data = my_data_tidy_proyecto)

mod2_c_soc <- lm(c_soc ~ 1 + conservatorio + sex, data = my_data_tidy_proyecto)

mod_int_c_soc <- lm(c_soc ~ 1 + conservatorio + sex + conservatorio:sex, data = my_data_tidy_proyecto)

summary(mod_int_c_soc)
## 
## Call:
## lm(formula = c_soc ~ 1 + conservatorio + sex + conservatorio:sex, 
##     data = my_data_tidy_proyecto)
## 
## Residuals:
##     Min      1Q  Median      3Q     Max 
## -5.6444 -1.4667  0.3556  1.3556  3.0000 
## 
## Coefficients:
##                         Estimate Std. Error t value Pr(>|t|)    
## (Intercept)               7.0000     0.2466  28.380   <2e-16 ***
## conservatoriosi           0.6444     0.3488   1.848   0.0664 .  
## sexmale                   0.6444     0.3488   1.848   0.0664 .  
## conservatoriosi:sexmale  -0.8222     0.4933  -1.667   0.0973 .  
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## 
## Residual standard error: 1.655 on 176 degrees of freedom
## Multiple R-squared:  0.0253, Adjusted R-squared:  0.008684 
## F-statistic: 1.523 on 3 and 176 DF,  p-value: 0.2103
anova(mod_null_c_soc, mod1_c_soc, mod2_c_soc, mod_int_c_soc)
## Analysis of Variance Table
## 
## Model 1: c_soc ~ 1
## Model 2: c_soc ~ 1 + conservatorio
## Model 3: c_soc ~ 1 + conservatorio + sex
## Model 4: c_soc ~ 1 + conservatorio + sex + conservatorio:sex
##   Res.Df    RSS Df Sum of Sq      F  Pr(>F)  
## 1    179 494.33                              
## 2    178 491.88  1    2.4500 0.8949 0.34544  
## 3    177 489.43  1    2.4500 0.8949 0.34544  
## 4    176 481.82  1    7.6056 2.7782 0.09734 .
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
plot(mod_int_c_soc)

durbinWatsonTest(mod_int_c_soc)
##  lag Autocorrelation D-W Statistic p-value
##    1      0.06295135      1.864671   0.316
##  Alternative hypothesis: rho != 0
mod_null_ef_fi <- lm(ef_fi ~ 1, data = my_data_tidy_proyecto)

mod1_ef_fi <- lm(ef_fi ~ 1 + conservatorio, data = my_data_tidy_proyecto)

mod2_ef_fi <- lm(ef_fi ~ 1 + conservatorio + sex, data = my_data_tidy_proyecto)

mod_int_ef_fi <- lm(ef_fi ~ 1 + conservatorio + sex + conservatorio:sex, data = my_data_tidy_proyecto)

summary(mod_int_ef_fi)
## 
## Call:
## lm(formula = ef_fi ~ 1 + conservatorio + sex + conservatorio:sex, 
##     data = my_data_tidy_proyecto)
## 
## Residuals:
##     Min      1Q  Median      3Q     Max 
## -3.8889 -0.9556  0.0444  1.0444  2.1111 
## 
## Coefficients:
##                         Estimate Std. Error t value Pr(>|t|)    
## (Intercept)              7.88889    0.19022  41.472   <2e-16 ***
## conservatoriosi          0.13333    0.26901   0.496    0.621    
## sexmale                  0.06667    0.26901   0.248    0.805    
## conservatoriosi:sexmale  0.20000    0.38044   0.526    0.600    
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## 
## Residual standard error: 1.276 on 176 degrees of freedom
## Multiple R-squared:  0.01427,    Adjusted R-squared:  -0.002528 
## F-statistic: 0.8496 on 3 and 176 DF,  p-value: 0.4685
anova(mod_null_ef_fi, mod1_ef_fi, mod2_ef_fi, mod_int_ef_fi)
## Analysis of Variance Table
## 
## Model 1: ef_fi ~ 1
## Model 2: ef_fi ~ 1 + conservatorio
## Model 3: ef_fi ~ 1 + conservatorio + sex
## Model 4: ef_fi ~ 1 + conservatorio + sex + conservatorio:sex
##   Res.Df    RSS Df Sum of Sq      F Pr(>F)
## 1    179 290.73                           
## 2    178 288.28  1      2.45 1.5047 0.2216
## 3    177 287.03  1      1.25 0.7677 0.3821
## 4    176 286.58  1      0.45 0.2764 0.5998
plot(mod_int_ef_fi)

durbinWatsonTest(mod_int_ef_fi)
##  lag Autocorrelation D-W Statistic p-value
##    1      0.07146747      1.838989   0.224
##  Alternative hypothesis: rho != 0
mod_null_leng_cas <- lm(leng_cas ~ 1, data = my_data_tidy_proyecto)

mod1_leng_cas <- lm(leng_cas ~ 1 + conservatorio, data = my_data_tidy_proyecto)

mod2_leng_cas <- lm(leng_cas ~ 1 + conservatorio + sex, data = my_data_tidy_proyecto)

mod_int_leng_cas <- lm(leng_cas ~ 1 + conservatorio + sex + conservatorio:sex, data = my_data_tidy_proyecto)

summary(mod_int_leng_cas)
## 
## Call:
## lm(formula = leng_cas ~ 1 + conservatorio + sex + conservatorio:sex, 
##     data = my_data_tidy_proyecto)
## 
## Residuals:
##     Min      1Q  Median      3Q     Max 
## -3.4222 -1.1111 -0.1111  0.9833  3.2667 
## 
## Coefficients:
##                         Estimate Std. Error t value Pr(>|t|)    
## (Intercept)               6.1111     0.2216  27.582  < 2e-16 ***
## conservatoriosi           1.3111     0.3133   4.184 4.51e-05 ***
## sexmale                   0.6222     0.3133   1.986   0.0486 *  
## conservatoriosi:sexmale  -0.6667     0.4431  -1.504   0.1343    
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## 
## Residual standard error: 1.486 on 176 degrees of freedom
## Multiple R-squared:  0.1175, Adjusted R-squared:  0.1025 
## F-statistic: 7.813 on 3 and 176 DF,  p-value: 6.321e-05
anova(mod_null_leng_cas, mod1_leng_cas, mod2_leng_cas, mod_int_leng_cas)
## Analysis of Variance Table
## 
## Model 1: leng_cas ~ 1
## Model 2: leng_cas ~ 1 + conservatorio
## Model 3: leng_cas ~ 1 + conservatorio + sex
## Model 4: leng_cas ~ 1 + conservatorio + sex + conservatorio:sex
##   Res.Df    RSS Df Sum of Sq       F    Pr(>F)    
## 1    179 440.58                                   
## 2    178 397.56  1    43.022 19.4751 1.773e-05 ***
## 3    177 393.80  1     3.756  1.7000    0.1940    
## 4    176 388.80  1     5.000  2.2634    0.1343    
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
plot(mod_int_leng_cas)

durbinWatsonTest(mod_int_leng_cas)
##  lag Autocorrelation D-W Statistic p-value
##    1      0.08570086      1.810181   0.178
##  Alternative hypothesis: rho != 0
mod_null_leng_cat <- lm(leng_cat ~ 1, data = my_data_tidy_proyecto)

mod1_leng_cat <- lm(leng_cat ~ 1 + conservatorio, data = my_data_tidy_proyecto)

mod2_leng_cat <- lm(leng_cat ~ 1 + conservatorio + sex, data = my_data_tidy_proyecto)

mod_int_leng_cat <- lm(leng_cat ~ 1 + conservatorio + sex + conservatorio:sex, data = my_data_tidy_proyecto)

summary(mod_int_leng_cat)
## 
## Call:
## lm(formula = leng_cat ~ 1 + conservatorio + sex + conservatorio:sex, 
##     data = my_data_tidy_proyecto)
## 
## Residuals:
##     Min      1Q  Median      3Q     Max 
## -3.3778 -1.3778  0.0222  1.0222  3.5556 
## 
## Coefficients:
##                         Estimate Std. Error t value Pr(>|t|)    
## (Intercept)               5.9778     0.2251  26.558  < 2e-16 ***
## conservatoriosi           1.4000     0.3183   4.398 1.89e-05 ***
## sexmale                   0.4667     0.3183   1.466    0.144    
## conservatoriosi:sexmale  -0.4667     0.4502  -1.037    0.301    
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## 
## Residual standard error: 1.51 on 176 degrees of freedom
## Multiple R-squared:  0.1415, Adjusted R-squared:  0.1269 
## F-statistic: 9.672 on 3 and 176 DF,  p-value: 6.08e-06
anova(mod_null_leng_cat, mod1_leng_cat, mod2_leng_cat, mod_int_leng_cat)
## Analysis of Variance Table
## 
## Model 1: leng_cat ~ 1
## Model 2: leng_cat ~ 1 + conservatorio
## Model 3: leng_cat ~ 1 + conservatorio + sex
## Model 4: leng_cat ~ 1 + conservatorio + sex + conservatorio:sex
##   Res.Df    RSS Df Sum of Sq       F    Pr(>F)    
## 1    179 467.39                                   
## 2    178 406.14  1     61.25 26.8664 5.935e-07 ***
## 3    177 403.69  1      2.45  1.0747    0.3013    
## 4    176 401.24  1      2.45  1.0747    0.3013    
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
plot(mod_int_leng_cat)

durbinWatsonTest(mod_int_leng_cat)
##  lag Autocorrelation D-W Statistic p-value
##    1      0.04466352      1.899911   0.438
##  Alternative hypothesis: rho != 0
mod_null_leng_est <- lm(leng_est ~ 1, data = my_data_tidy_proyecto)

mod1_leng_est <- lm(leng_est ~ 1 + conservatorio, data = my_data_tidy_proyecto)

mod2_leng_est <- lm(leng_est ~ 1 + conservatorio + sex, data = my_data_tidy_proyecto)

mod_int_leng_est <- lm(leng_est ~ 1 + conservatorio + sex + conservatorio:sex, data = my_data_tidy_proyecto)

summary(mod_int_leng_est)
## 
## Call:
## lm(formula = leng_est ~ 1 + conservatorio + sex + conservatorio:sex, 
##     data = my_data_tidy_proyecto)
## 
## Residuals:
##     Min      1Q  Median      3Q     Max 
## -4.0667 -1.6000  0.2444  1.2667  2.9333 
## 
## Coefficients:
##                         Estimate Std. Error t value Pr(>|t|)    
## (Intercept)               7.0667     0.2558  27.627   <2e-16 ***
## conservatoriosi           0.7111     0.3617   1.966   0.0509 .  
## sexmale                   0.5333     0.3617   1.474   0.1422    
## conservatoriosi:sexmale  -0.5778     0.5116  -1.129   0.2603    
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## 
## Residual standard error: 1.716 on 176 degrees of freedom
## Multiple R-squared:  0.02716,    Adjusted R-squared:  0.01058 
## F-statistic: 1.638 on 3 and 176 DF,  p-value: 0.1823
anova(mod_null_leng_est, mod1_leng_est, mod2_leng_est, mod_int_leng_est)
## Analysis of Variance Table
## 
## Model 1: leng_est ~ 1
## Model 2: leng_est ~ 1 + conservatorio
## Model 3: leng_est ~ 1 + conservatorio + sex
## Model 4: leng_est ~ 1 + conservatorio + sex + conservatorio:sex
##   Res.Df    RSS Df Sum of Sq      F Pr(>F)
## 1    179 532.64                           
## 2    178 524.62  1    8.0222 2.7248 0.1006
## 3    177 521.93  1    2.6889 0.9133 0.3406
## 4    176 518.18  1    3.7556 1.2756 0.2603
plot(mod_int_leng_est)

durbinWatsonTest(mod_int_leng_est)
##  lag Autocorrelation D-W Statistic p-value
##    1      0.03697763      1.916163   0.524
##  Alternative hypothesis: rho != 0
mod_null_mat <- lm(mat ~ 1, data = my_data_tidy_proyecto)

mod1_mat <- lm(mat ~ 1 + conservatorio, data = my_data_tidy_proyecto)

mod2_mat <- lm(mat ~ 1 + conservatorio + sex, data = my_data_tidy_proyecto)

mod_int_mat <- lm(mat ~ 1 + conservatorio + sex + conservatorio:sex, data = my_data_tidy_proyecto)

summary(mod_int_mat)
## 
## Call:
## lm(formula = mat ~ 1 + conservatorio + sex + conservatorio:sex, 
##     data = my_data_tidy_proyecto)
## 
## Residuals:
##     Min      1Q  Median      3Q     Max 
## -4.5333 -1.2444  0.0889  1.3778  3.7556 
## 
## Coefficients:
##                         Estimate Std. Error t value Pr(>|t|)    
## (Intercept)               5.5333     0.2779  19.910  < 2e-16 ***
## conservatoriosi           2.0889     0.3930   5.315  3.2e-07 ***
## sexmale                   0.7111     0.3930   1.809   0.0721 .  
## conservatoriosi:sexmale  -0.4222     0.5558  -0.760   0.4485    
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## 
## Residual standard error: 1.864 on 176 degrees of freedom
## Multiple R-squared:  0.2194, Adjusted R-squared:  0.2061 
## F-statistic: 16.49 on 3 and 176 DF,  p-value: 1.737e-09
anova(mod_null_mat, mod1_mat, mod2_mat, mod_int_mat)
## Analysis of Variance Table
## 
## Model 1: mat ~ 1
## Model 2: mat ~ 1 + conservatorio
## Model 3: mat ~ 1 + conservatorio + sex
## Model 4: mat ~ 1 + conservatorio + sex + conservatorio:sex
##   Res.Df    RSS Df Sum of Sq       F   Pr(>F)    
## 1    179 783.66                                  
## 2    178 624.99  1   158.672 45.6511 1.99e-10 ***
## 3    177 613.74  1    11.250  3.2367  0.07372 .  
## 4    176 611.73  1     2.006  0.5770  0.44850    
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
plot(mod_int_mat)

durbinWatsonTest(mod_int_mat)
##  lag Autocorrelation D-W Statistic p-value
##    1       0.1154332      1.734837   0.068
##  Alternative hypothesis: rho != 0
mod_null_plast <- lm(plast ~ 1, data = my_data_tidy_proyecto)

mod1_plast <- lm(plast ~ 1 + conservatorio, data = my_data_tidy_proyecto)

mod2_plast <- lm(plast ~ 1 + conservatorio + sex, data = my_data_tidy_proyecto)

mod_int_plast <- lm(plast ~ 1 + conservatorio + sex + conservatorio:sex, data = my_data_tidy_proyecto)

summary(mod_int_plast)
## 
## Call:
## lm(formula = plast ~ 1 + conservatorio + sex + conservatorio:sex, 
##     data = my_data_tidy_proyecto)
## 
## Residuals:
##     Min      1Q  Median      3Q     Max 
## -3.6000 -0.9556 -0.1778  0.8222  2.4000 
## 
## Coefficients:
##                         Estimate Std. Error t value Pr(>|t|)    
## (Intercept)               6.9556     0.1982  35.097  < 2e-16 ***
## conservatoriosi           1.2222     0.2803   4.361  2.2e-05 ***
## sexmale                   0.6444     0.2803   2.299   0.0227 *  
## conservatoriosi:sexmale  -0.6000     0.3964  -1.514   0.1319    
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## 
## Residual standard error: 1.329 on 176 degrees of freedom
## Multiple R-squared:  0.1329, Adjusted R-squared:  0.1181 
## F-statistic: 8.989 on 3 and 176 DF,  p-value: 1.429e-05
anova(mod_null_plast, mod1_plast, mod2_plast, mod_int_plast)
## Analysis of Variance Table
## 
## Model 1: plast ~ 1
## Model 2: plast ~ 1 + conservatorio
## Model 3: plast ~ 1 + conservatorio + sex
## Model 4: plast ~ 1 + conservatorio + sex + conservatorio:sex
##   Res.Df    RSS Df Sum of Sq       F  Pr(>F)    
## 1    179 358.73                                 
## 2    178 320.46  1    38.272 21.6542 6.4e-06 ***
## 3    177 315.12  1     5.339  3.0207 0.08396 .  
## 4    176 311.07  1     4.050  2.2915 0.13188    
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
plot(mod_int_plast)

durbinWatsonTest(mod_int_plast)
##  lag Autocorrelation D-W Statistic p-value
##    1     -0.09708371      2.171191    0.32
##  Alternative hypothesis: rho != 0
mod_null_tecn <- lm(tecn ~ 1, data = my_data_tidy_proyecto)

mod1_tecn <- lm(tecn ~ 1 + conservatorio, data = my_data_tidy_proyecto)

mod2_tecn <- lm(tecn ~ 1 + conservatorio + sex, data = my_data_tidy_proyecto)

mod_int_tecn <- lm(tecn ~ 1 + conservatorio + sex + conservatorio:sex, data = my_data_tidy_proyecto)

summary(mod_int_tecn)
## 
## Call:
## lm(formula = tecn ~ 1 + conservatorio + sex + conservatorio:sex, 
##     data = my_data_tidy_proyecto)
## 
## Residuals:
##     Min      1Q  Median      3Q     Max 
## -3.5778 -1.0667 -0.0667  0.9333  2.7778 
## 
## Coefficients:
##                         Estimate Std. Error t value Pr(>|t|)    
## (Intercept)               6.5778     0.2067  31.824  < 2e-16 ***
## conservatoriosi           1.4889     0.2923   5.094 8.99e-07 ***
## sexmale                   0.6444     0.2923   2.205   0.0288 *  
## conservatoriosi:sexmale  -0.6444     0.4134  -1.559   0.1208    
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## 
## Residual standard error: 1.387 on 176 degrees of freedom
## Multiple R-squared:  0.1726, Adjusted R-squared:  0.1585 
## F-statistic: 12.24 on 3 and 176 DF,  p-value: 2.595e-07
anova(mod_null_tecn, mod1_tecn, mod2_tecn, mod_int_tecn)
## Analysis of Variance Table
## 
## Model 1: tecn ~ 1
## Model 2: tecn ~ 1 + conservatorio
## Model 3: tecn ~ 1 + conservatorio + sex
## Model 4: tecn ~ 1 + conservatorio + sex + conservatorio:sex
##   Res.Df    RSS Df Sum of Sq       F    Pr(>F)    
## 1    179 408.95                                   
## 2    178 347.70  1    61.250 31.8600 6.524e-08 ***
## 3    177 343.03  1     4.672  2.4303    0.1208    
## 4    176 338.36  1     4.672  2.4303    0.1208    
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
plot(mod_int_tecn)

durbinWatsonTest(mod_int_tecn)
##  lag Autocorrelation D-W Statistic p-value
##    1       0.2065035       1.57429   0.004
##  Alternative hypothesis: rho != 0